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We investigate the properties of triply ionized Carbon (CIV) in the Intergalactic Medium 
(IGM) using a set of high-resolution and large box-size cosmological hydrodynamic simula- 
tions of a ACDM model. We rely on a modification of the public available Tree-Smoothed 
Particle Hydrodynamics (SPH) code GADGET-2 that self-consistently follows the metal en- 
richment mechanism by means of a detailed chemical evolution model. We focus on several 
numerical implementations of galactic feedback: galactic winds in the energy driven and mo- 
mentum driven prescriptions, galactic winds hydrodynamically coupled to the surrounding 
gas and Active Galactic Nuclei (AGN) powered by gas accretion onto massive black holes. 
Furthermore, our results are compared to a run in which galactic feedback is not present and 
we also explore different initial stellar mass function. After having addressed some global 
properties of the simulated volume like the impact on the star formation rate and the content in 
Carbon and CIV, we extract mock IGM transmission spectra in neutral hydrogen (HI) and CIV 
and perform Voigt profile fitting. The results are then compared with high-resolution quasar 
(QSO) spectra obtained with the UVES spectrograph (Ultra Violet Echelle Spectrograph) at 
the VLT (Very Large Telescope) and the High Resolution Echelle Spectrometer (HIRES) at 
Keck. We find that feedback has little impact on statistics related to the neutral hydrogen, 
while CIV is more affected by galactic winds and/or AGN feedback. The feedback schemes 
investigated have a different strength and redshift evolution with a general tendency for AGN 
feedback to be more effective at lower redshift than galactic winds. When the same analysis 
is performed over observed and simulated CIV lines, we find reasonables good agreement 
between data and simulations over the column density range A/civ = 10 125 ~ 15 cm~ 2 . Also 
the CIV line-widths distribution appears to be in agreement with the observed values, while 
the HI Doppler parameters, 6hi, are in general too large showing that the diffuse cosmic web 
is heated more than what is inferred by observations. The simulation without feedback fails 
in reproducing the CIV systems at high column densities at all redshift, while the AGN feed- 
back case agrees with observations only at z < 3, when this form of feedback is particularly 
effective. We also present scatter plots in the b — N and in the New — A^hi planes, showing 
that there is rough agreement between observations and simulations only when feedback is 
taken into account. 

Although it seems difficult to constrain the nature and the strength of galactic feedback 
using the present framework and find a unique model that fits all the observations, these sim- 
ulations offer the perspective of understanding the galaxy-IGM interplay and how metals pro- 
duced by stars can reach the low density IGM. 

Key words: methods: numerical - galaxies: formation - intergalactic medium - quasars: 
absorption lines - cosmology: theory 
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1 INTRODUCTION 

The interplay between galaxies and the intergalactic medium 
(IGM) is a fundamental and poorly understood aspect in the stan- 
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dard scenario of structure formation. The intimate connection be- 
tween the diffuse IGM and the baryons that populate star forming 
galaxies can shed light on the metal enrichment of the IGM (its na- 
ture and redshift evolution), on the properties of galactic feedback 
and on the evolution of the ultraviolet background. 

In recent years, this field has benefitted from spectroscopi- 
cal observations of quasar (QSO) and Gamma Ray Burst (GRB) 
spectra (both at high and low resolution) and from photometric 
and spectroscopical information on high redshift galaxies in the 
same field of view. This has in turn allowed us to perform cross- 
correl ation studies of a bsorption lines and properties of galaxies 
(e.g. lSteidel et alj20ld . and references therein). Observations have 
revealed that the IGM is pollut ed with heavy elements such as Car- 
bon, Oxygen and Silicon (e.g. IWevmann et alj|l98ll ; ICowie et al.l 
ll995[|Aguirre et alj2005l) . Furthermore, high redshift star forming 
galaxies have been shown to interact with t he surrounding medium 
by powering strong galactic outflows (e.g. lAdelberger et alj|2005l) . 
The IGM enrichment mechanism is, however, still poorly con- 
strained and the extent of metal pollution of the low d ensity IGM, 
or around Lyman-break galaxies, controversial (e.g. ISchave et al.l 
l2000al : |PierietalJ [2006: Por ciani & Madaull2005l) . 

Among many different metal species, the doublet produced 
by triply ionized carbon (CIV) is the one that has been most widely 
studied since it is relatively easy to identify in the regions redward 
of Lyman-a emission. CIV has been obse rved from low-redshift up 
to ver y high redshift flRvan- Weber et alj|20 06. 2009; Becker et al.l 
|2009J); its clustering properties in the transv erse direction and along 
the line-of-sight have been investiga ted JD'Odorico et al.l 120021 ; 
iPichon et al]|2003l : iTvtler et alj|2009t) along wit h the global evo- 
lution of the CIV density fraction down to z < 2 dD'Odorico et al.l 
l20ld) ; the optical depth in CIV has been derived using pixel- 
optical-depth techniques applied to the UVES Large Programme 
sample in order to constrain the metallicity of the low density gas 
dSchave et a l. 2003; Aracil et al. 2004); the cross-correlation be- 
tween galaxies and CIV properties have also been addressed in 
great detail using either K eck or UVES spectra dSimcoe et al . 2006; 
IScannapieco et alj|2006l) . The general conclusion is that CIV ab- 
sorption systems are indeed related to outskirts of galactic haloes 
but it is still unclear whether the gas producing the absorption is en- 
riched by nearby galaxies or at an earlier epoch by smaller objects. 

The large amount of data available has prompted the imple- 
mentation of relevant physical processes into sophisticated hydro- 
dynamic codes in order to simulate both galaxies and QSO/GRB 
absorption lines. The first attempt to model CIV properties us- 
ing hy drodynamic simulations has been made by lHaehnelt et al.l 
( 1996). Subsequently, it has been demonstrated that a key in- 
gredient that appears to be necessary in order to reproduce CIV 
metal line statistics is feedback in the form of galactic winds 
dTheuns et alj|2002l ; ICen et alj|2005l ; lOppenheimer & Davell2006l 
2008). Galactic winds triggered by supernovae (SN) explosions 
have been modelled in a simple way using mainly two schemes that 
conserve either ene rgy or momentum (Oppenheimer & Dave 2006; 
iTescari et al.l2009l) . Other sources of feedback like Active Galactic 
Nuclei (AGN) has also been cons idered to address glob al gas prop- 
erties from high to low redshift l lTornatore et alJuOlOl) . A second 
very important ingredient has been the implementation of accurate 
chemical evolution models that self-consistently follow the release 
of different metal species fr om stars and also allow for different ini- 
tial st e llar mass function s 1 Tornatore et alj|2007l; Kawata & Rauch 
20071; iKobavashi et alj |2007|; ISommer-Larsen & Fvnbol l200a 



Wiersma et al. 2009; Oppenheimer et al. 2009; Wiersmaetal. 
20ld;lDave et al.ll2010l ; lcen & Chisarill2010h . 



This is a second paper in a series that exploits the capabili- 
ties of a state-of-the-art modification of the public available code 
GADGET-2 (Springel 2005). In the first paper we addressed neutral 
hydrogen properties JTescari et al.ll2009l) focusing in particular on 
Damped Lyman-a systems. In this paper we study the CIV transi- 
tion by extending the feedback models simulated. 

This paper is organized as follows. In Sections [2] and [3] we 
present, respectively, the observational data sample and our set of 
simulations along with the different feedback models implemented 
in our code. In Section|4]we compare some global properties of the 
simulations: the cosmic star formation rate and the evolution of the 
total Carbon content. We then focus on the neutral hydrogen statis- 
tics in Section[5] In Sections[6]and|7]we investigate the CIV column 
density distribution function and the evolution with redshift of the 
CIV cosmological mass density, f^civ, respectively. In Section [8] 
we study the probability distribution function of the CIV Doppler 
parameter and, in Section [9] we present the CIV column density- 
Doppler parameter relation. The last part of the paper is dedicated 
to the analysis of the HI-CIV correlated absorption (Section [Tot. 
Finally, in Section[TT|we summarize our main results and we draw 
some conclusions. 



2 OBSERVATIONAL DATA SAMPLE 

The core of the observational data sample is formed by the high 
resolution (R ~ 45000), high signal-t o-noise (S/N) QSO spectra 
obtained with the UVES spectrograph JDekker et alJ uOOO) at the 
Kueyen unit of the ESO VLT (Cerro Paranal, Chile) in the frame- 
work of the ESO Large Programm e (LP): "The Cosmic Evolution 
of the IGM" dBergeron et alj2004) . 

We ha ve used the sampl e of HI Lyman-a absorption l ines de- 
scribed in lSaittaetaT] <2008h and iD'Odorico et alj J2008I) . which 
is based on 18 objects of the LP plus 4 more QSOs of which: 
J2233-606 and HE1 122-1648 were observed with UVES, while 
HS 1946+7658 and B 1422+231 were obtained with the High Res- 
olution Echelle Spectrometer (HIRES) at the Keck Telescope at 
comparable resolution and signal-to-noise ratio. In this work, also 
the HI Lyman-a list of the QSO Q0055-269, observed with UVES, 
has been added to the sample. The adopted s ample of CIV absorp- 
tions was described in lD'Odorico et alj J20101) . The previous group 
of QSOs was increased with the addition of the UVES spectra of 
PKS2000-330 and PKS 1937-101 to extend the redshift interval to 

Z ru O. 

All UVES spectra were reduced with the UVES pipeline fol- 
lowing the standard procedure. The continuum level was deter- 
mined by interpolating with a cubic spline the region of the spec- 
trum free from evident absorption features. Absorption lines were 
fitted with Voigt pr ofiles using the conte xt LYMAN of the MIDAS 
reduction package JFontana & Ballesterlll995l) . A minimum num- 
ber of components was adopted to fit the velocity profile in order 
to reach a normalized \ 2 ~ 1- F° r a U me QSO in the sample, the 
Lyman-a forest was defined as the interval between the Lyman- 
/3 emission plus 1000 km s _1 (to avoid contamination from the 
Lyman-/3 forest), and 5000 km s _1 from the Lyman-a emission; 
the CIV forest was defined as the interval between the Lyman- 
a emission and 5000 km s _1 from the CIV emission. The higher 
bound was fixed to avoid the proximity region affected by the QSO 
emission, where most of the intrinsic systems are found. 

CIV features are distinguished into components or lines, 
meaning the velocity components in which every absorption pro- 
file has been decomposed, and systems, formed by groups of com- 
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QSO 


^cm 


AzLya 


Azciv 


HE 1341-1020 


2.142 


1.6599-2.090 


1.467-2.090 


Q0122-380 
PKS 1448-232 


2.2004 

2.224 


1.709-2.147 
1.729-2.171 


1.513-2.147 
1.531-2.171 


PKS 0237-23 


2.233 


1.737-2.179 


1.538-2.179 


J2233-606 
HE 0001-2340 


2.248 
2.265 


1.7496-2.194 
1.764-2.211 


1.550-2.194 
1.564-2.211 


HS 1626+6433" 
HE1122-1648 6 


2.32 
2.40 


1.878-2.344 


1.607-2.265 
1.665-2.344 


Q0 109-35 18 
HE 2217-2818 


2.4057 
2.414 


1.883-2.349 
1.888-2.355 


1.674-2.349 
1.681-2.355 


Q0329-385 
HE 1158-1843 


2.435 
2.448 


1.9096-2.378 
1.919-2.391 


1.697-2.378 
1.707-2.391 


HE 1347-2457 
Q 1442+2931" 
Q0453-423 
PKS 0329-255 


2.5986 
2.661 
2.669 
2.696 


2.046-2.539 

2.106-2.608 
2.129-2.635 


1.826-2.539 
1.875-2.600 
1.881-2.608 
1.902-2.635 


HE 0151-4326 
Q0002-422 
HE 2347-4342 
SBS 1107+487" 


2.763 
2.769 
2.880 
2.966 


2.186-2.701 
2.191-2.707 
2.285-2.816 


1.955-2.701 
1.959-2.707 
2.067-2.816 
2.114-2.900 


HS 1946+7658 c 
HE 0940-1050 


3.058 
3.0932 


2.435-2.991 
2.465-3.025 


2.181-2.991 
2.214-3.025 


QO420-388 
S4 0636+68" 


3.1257 
3.175 


2.493-3.057 


2.239-3.057 
2.278-3.106 


SBS 1425+606" 
PKS 2126-158 


3.199 
3.292 


2.633-3.221 


2.297-3.129 
2.370-3.221 


B 1422+23 l b 
Q0055-269 
PKS 2000-330 
PKS 1937-101 


3.623 
3.66 
3.783 
3.787 


2.914-3.546 
2.945-3.583 


2.630-3.546 
2.659-3.583 
2.756-3.704 
2.770-3.400 


PSS J 1646+55 14" 


4.059 




2.972-3.975 


PSS J1057+4555" 
BR 2237-0607" 


4.131 
4.559 




3.029-4.046 
3.365-4.467 
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pon ents with relative separation smaller than du m i n = 50 km s _1 
(see ID'Odorico et alJuOld . for a detailed description of the proce- 
dure). Both the velocity components and the systems in the sample 
have column densities in the range 10 11 ^ iVciv ~ 10 15 cm~ 2 . 
In order to further increase the number of CIV lines and to 

• 1 j — 7y 

D'Odorico et al. (2010) have 
considered the CIV lines fitted in 9 QSO spectra observed with 
HIRES at Keck at a resolution and si gnal-to-noise ratio simila r 
to those of their spectra and reported in Boksenberg et al. (2003). 
The fi t with Voigt profiles was carried out by Boksenberg et al. 
(2003) with VPFnQ. The main difference between LYMAN and 
VPFIT is that the number of components fitted to a given veloc- 
ity profile i s, in general, lar ger using the latter (see also the dis- 
cussion in Saitt a et alj 120080 . This is seen also in the analysis 
of ID'Odorico et alj J2010h . in particular from the comparison of 
the CIV lines detected in the spectrum of the QSO B 1422+231, 
which is the on l y obj ect in common between the two samples. 
ID'Odorico et alj 12010]) find that in all cases the number of compo- 
nents determined with VPFIT is larger or equal to that determined 
with LYMAN. However, when the total column density of each ab- 
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sorption system is considered the difference between the two fitting 
procedures becomes negligible. 

All the QSOs of the observational sample are reported in Table 
1 with their emission redshift and the redshift range covered by the 
Lyman-Q forest and the CIV lines. 



3 THE SIMULATIONS 

In the following we review the mai n characteristi cs of the runs anal- 
ysed in this work and we refer to iTescari et alj 120090 for a more 
extensive description of the simulations. Our code i s the same mod- 
ified v ersion of the TreePM -SPH code GADGET-2 JSpringelll20o3) 
used in ITescari et alj 12009J). 

The simulations cover a cosmological volume (with periodic 
boundary conditions) filled with an equal number of dark matter 
and gas particles. The cosmological model chosen is a flat ACDM 
with the following parameters: ilo m = 0.24, Slob = 0.0413, 
tt A = 0.76, n s = 0.96, Ho = 73 km s" 1 Mpc" 1 and <t 8 = 0.8, 
which are in agreement with the latest results from large scale struc- 
ture observables such as the cosmic microwave background, weak 
lensing, the Lyma n-Q forest and the evolu t ion of mass function 
of galaxy clusters JLesgourgues et alj|2007l ; iKomatsu et alj|2009t 
IVikhlinin et al. 2009). The input linear dark matter power spectrum 
for the initial conditions has been generated at z — 99 with CMB- 
FAST (Seli ak & Zaldarriaga 1996) and includes baryonic acoustic 
oscillations. 

Radiative cooling and heating processes are modelled follow- 
ing the implementation of Katz et al. ( 1996). We assume a mean 
Ultra Violet Background (UVB) p roduced by quasars and galaxies 
as given by IHaardt & Madaul (1996J), with the heating rates multi- 
plied by a factor of 3.3 in order to better fit observational constraints 
on the temperature evolution of the intergalactic medium at high 
redshift (the factor of 3.3 was originally introduced because it gives 
a temperature at the mean cosmic density, To, in agreement with the 
ISchave et alj|2000bl . data points). This background gives naturally 
r_i2 ~ 0.8 — 1 at lo w redshift, in agreement with recent obser- 
vational measurements jBolton et al.ll2005l ; iFaucher-Giguere et alj 
2008). Multiplying the heating rates by the factor above (chosen 
empirically) results in a larger IGM temperature at the mean den- 
sity which cannot be reached by the standard hydrodynamic code. 
This aims at mimicking, at least in a phenomenological way, the 
non-equilibrium ionization effects around reionization (see for ex- 
ample iRicotti et alj[200ol: ISchave et ai]|2000al : iBolton & Haehneld 
120071) . A significant fraction in mass of the IGM resides in the tight 
power-law relation T = Tb(l + <5) 7_1 at around the mean den- 
sity, where To is the temperature at the mean cosmic density and 
5 the overdensity. Our reference simulations (see at the end of this 
Section), at redshift z — 3, 2.25 and 1.5, have log{Tb(K)}, respec- 
tively, equal to: 4.36, 4.33, 4.30 and 7, respectively, equal to: 1.58, 
1.60, 1.62. 

In our simulations, the standard m ultiphase star formation cri- 
terion of lSpringel & Hernquistl (2003J) is implemented, which con- 
tains an effective prescription for the inter-stellar medium (ISM). In 
this effective model, whenever the density of a gas particle exceeds 
a given threshold (set to pth = 0.1 cm -3 in terms of the num- 
ber density of hydrogen atoms), that gas particle is flagged as star 
forming and is treated as multiphase. With this prescription baryons 
are in the form either a hot or a cold phase or in stars, thereby this 
density threshold marks the onset of cold clouds formation. 

We follow self-consistently the evolution of these elements: H, 
He, C, O, Mg, S, Si and Fe. The contribution of metals is included 
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Run 


IMF 


Box Size 
[/i _1 Mpc] 


Ngas 


mGAS 


Softening 
[/i -1 kpc] 


Feedback 


kr37edw500 


Kroupa 


37.5 


256 3 


3.6 x 10 7 


7.5 


EDW 


kr37mdw 


Kroupa 


37.5 


256 3 


3.6 x 10 7 


7.5 


MDW 


kr37agn 


Kroupa 


37.5 


256 3 


3.6 x 10 7 


7.5 


AGN 6 


kr37agn+edw300 


Kroupa 


37.5 


256 3 


3.6 x 10 7 


7.5 


AGN+EDW C 


kr37nf 


Kroupa 


37.5 


256 3 


3.6 x 10 7 


7.5 


NO FEEDBACK 


kr37co-edw500 


Kroupa 


37.5 


256 3 


3.6 x 10 7 


7.5 


COUPLED EDW d 


ay37edw500 


Arimoto-Yoshii 


37.5 


256 3 


3.6 x 10 7 


7.5 


EDW 


sa37edw500 


Salpeter 


37.5 


256 3 


3.6 x 10 7 


7.5 


EDW 


kr37p400edw500 


Kroupa 


37.5 


400 3 


9.4 x 10 6 


4.8 


EDW 



Table 2. Summary of the different runs. Column 1, run name; column 2, Initial Mass Function (IMF) chosen; column 3, comoving box size; column 4, number 
of gas particles; column 5, mass of gas particle; column 6, Plummer-equivalent comoving gravitational softening; column 7, type of feedback implemented. For 
the runs kr37edw500, kr37co-edw500, ay37edw500, sa37edw500 and kr37p400edw500 the velocity of the energy-driven winds (EDW) is set to v w = 500 km 
s _1 . (a): momentum-driven winds (MDW) feedback in which wind velocity scales roughly with 3a (a being the velocity dispersion of the halo that hosts the 
"wind" particle), (b): AGN feedback in which the energy is released by gas accretion onto super-massive black holes, (c): combined effect of energy-driven 
winds of 300 km s~ 1 and AGN feedbacks, (d): energy-driven winds not decoupled from the hydrodynamics (see Section l3~lY 



in the cooling function using the tables of ISutherland & Dopital 
( 1993), that assume the solar value for the relative abundances. 
In this paper w e use the sola r meta llicity and element abun- 
dances given by lAsplund et all J2005D . with Z/X = 0.0165. 
Besides including different contributions from Type II and Type 
la supernovae (SN II, SN la) and Low and Intermediate-Mass 
Stars (LIMS), our model of chemical evolution accounts for mass- 
dependent stellar lifetimes. W e adopt the lifetime function given by 
Padovani & Matteucci 1 19931) . F or the stellar yields we use: SNIa - 



Thielemann et al.N2003l); SNII -IWooslev & Weaver! J 1995b ; LIMS 
- Ivan d en Hoek & Gr oenewegenl ( 119970 . The mass-range for the 
SNII ism > 8M Q , while for the SNIa is m < 8M Q with a 
binary fraction of 10%. Finally we use three distinct stellar ini- 
tial mass functions (IMFs): a Kroupa, a Salpeter and an Arimoto- 
Yoshi i IMF. For t his paper our reference choice is the Kroupa 
IMF d Kroup a 2001), which consists of a multi-slope approxima- 
tion: j/kr = 0.3 for stellar mass m < 0.5Mq, j/kr = 1.2 for 
0.5M Q < m < 1M Q and j/ K r = 1.7 for m > 1M Q (where 
f{rn) oc m _y is the functi onal form for the IMF and <p(m) is 
the IMF by mass). Salpeter (Salpeter 1955) IMF has single slope 
j/sl = 1-35 and Arimoto-Yoshii {Arimoto & Yoshii 1987) IMF 
also has single slope j/ay = 0.95. 

In this paper we use three different feedback schemes: two 
associated to the galactic (energy and momentum driven) winds 
produ ced by "starburst" galaxies, already presented in lTescari et al.l 
(2009), and the (new) AGN feedback associated to the energy re- 
leased by gas accretion onto super-massive black holes. In the fol- 
lowing we briefly describe these different models. 



Star forming gas particles are then stochastically selected accord- 
ing to their star formation rate to become part of a blowing wind. 
Whenever a particle is uploaded to the wind, it is decoupled from 
the hydrodynamics for a given period of time or till the density 
around it drops below a given fraction of the density threshold for 
the onset of the star formation (p t h = 0.1 cm -3 in terms of the 
number density of hydrogen atoms), in order to effectively reach 
less dense regions. The (lower) limiting density for decoupling is 
set to 0.5pth- This allows the wind particle to travel "freely" up 
to few kpc until it has left the dense star forming phase, without 
directly interacting with it. Otherwise, to verify the effect of the de- 
coupling, we run a "coupled winds" simulation in which wind parti- 
cles are always aff ected by the hydrodynamics like all the other gas 
partic les (see also Dalla Vecchia & Schave 2008]; iTornatore et al.l 
uOld) . Unlike in Springel & Hernquisn ( 120031) ■ we fixed the veloc- 
ity of the winds, w w , instead of the fraction of the energy made 
available by SNII explosions to power galactic ejecta. Thus, four 
parameters fully specify the wind model: the wind efficiency r\, the 
wind speed v w , the wind free travel length Z w and the wind free 
travel density factor 5„. The maximum allowed time for a wind 
particle to stay hydrodynamically decoupled is then tdcc = i w /"w 
(note that in the case of coupled winds £d oc = 0)- The parameter 
Z w has been introduced in order to prevent a gas particle from being 
trapped into the potential well of the virialized halo and in order 
to effectively escape from the ISM, reach the low density IGM and 
pollute it with metals. 

In this paper, we use Z w = 10 ft" 1 kpc and we considered two 
different values for the wind velocity: u w = 300 and 500 km s _1 . 
In our implementation the parameter r\ is kept fixed to the value 2. 



3.1 Energy-driven winds 

For this implementation, the wind mass-loss rate M w is assumed 
to be proportional to the star formation rate M* according to 



M w = rjM* 



(1) 



where 77 is the wind mass loading factor (i.e. the wind efficiency), 
and the wind carries a fixed fraction y of SN energy: 



1-2 

-M w u w = ycsnM* 



(2) 



3.2 Momentum-driven winds 

A feasib le physical sc e nario for the wind driving mechanism, out- 
lined by Murray et al. (2005) to explain observational evidence of 
outflows from starburst galaxies (Martin 2005; Rupke et al. 2005), 
is the so called momentum-driven wind model. In this scenario, 
the radiation pressure of the starburst drives an outflow, possibly 
by transferring momentum to an absorptive component (e.g. dust) 
that is hydrodynamically coupled with the gas compo nent, which 
is then dragged out of the galaxy. Following Opp enheimer & Pavel 
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(2006, 2008) we implemented the momentum-driven winds model 
in our code. In such a model the wind speed scales as the galaxy 
velocity dispersion a, as observed by Martin (2005). Since in 
momentum-driven winds the amount of input momentum per unit 
star formation is constant, this implies that the wind mass loading 
factor r\ must be inversely proportional to the velocity dispersion. 
We therefore use the following relations: 



^wind 



■Sa-JK 



V= — , 
a 



(3) 



where /l is the luminosity factor in units of the galactic Eddington 
luminosity (i.e. the critical luminosity necessary to expel gas from 
the galaxy), <jq is the normalization of the mass loading factor and 
we add an extra 2a kick to get out of the potential of the galaxy 
in order to simulate continuous pumping of gas until it is well into 
the galactic halo. To determine a, we identify haloes runtime in the 
simulations by using a parallel friends-of-friends (FoF) algorithm 
and we associate a with the velocity dispersion of the ha l oes, udm 
(for fu rther details see Section 2.2 of lTescari et alJl2009l) . ISchaererl 
1 2003) suggeste d an approximate functional for m for far-UV emis- 
sion which lOppenheimer & Dave (2006, 2008) used to obtain the 
luminosity factor and which includes a metallicity dependence for 
/l, owing to more UV photons output by lower- metallicity stellar 
populations: 



-0.0029(log Z+9) zo +0. 417694 



(4) 



/l = ,/l,o x 10 

where Martin (2005) suggests: /l,q ~ 2. The mass loading fac- 
tor controls star formation at early times, so gq can also be set by 
requiring a matc h to the obse r ved g l obal s tar formation rate. Fol- 
lowing Op penheimer & Pavel J2006I. [2008) we set cr = 150 km 
s -1 . Even for this wind implementation the particles are stochasti- 
cally selected in the same way as for the energy-driven scenario. 



effici ent accretion onto a Schwarzschild BH IShakura & Sunvaevl 
119730 , The model then assumes that a fraction et of the radiated 
energy is thermally coupled to the surrounding gas, so that 
E{ cc d = e r efA^BHC 2 is the rate of the e nergy released t o hea t 
the surrounding gas. Using et ~ 0.05, |Pi Matteo et alj (2005) 
were able to reproduce the observed Mbh — & relation between 
bulge velocity disper sion and mass of t he hosted BH (see also 
ISiiacki et alll2008l: |Pi Matteo et alj|2008h . Gas particle accretion 
onto the BH is implemented in a stochastic way, by assigning 
to each neighbouring gas particle a probability of contributing 
to the accretion, which is proportional to the SPH kernel weight 
computed at the particle position. In the scheme described above, 
this stochastic accretion is used only to increase the dynamic mass 
of the BHs, while the mass entering in the computation of the 
accretion rate is followed in a continuous way, by integrating the 
analytic expression for Mbh- Once the amount of energy to be 
thermalized is computed for each BH at a given timestep, this 
energy is then distributed to the surrounding gas particles using 
the SPH kernel weighting scheme. Our research group explored 
for the first time this new feedback prescription in relation to 
the high redsh i ft prop erties of the intergalactic medium, while 
iTornatore et al.1 ( 120100 have studied the impact of this feedback 
mechanism on the low redshift IGM (other groups are working 
on the black holes accretion feedback related to IGM and the 
Warm Hot Intergala ctic Medium (W HIM), see for example 
IWiersma eta l. 2010; Bertone et alj|201(fl ). 

In Table [2] we summarize the main parameters of the cosmo- 
logical simulations performed, including the mass associated to the 
gas particles and the gravitational softening. All the simulations 
start at redshift z — 99 and stop at redshift z = 1.5. In the follow- 
ing we outline the specific characteristics of each run: 



3.3 AGN feedback 

Finally, we include in our simulations the effect of feedback energy 
from gas accretion onto super-massive black holes (BHs), f ollow- 
ing th e scheme originally i ntroduced bv lSpringel et al.1 J20050 (see 
also lDi Matteo et alj|2005h . In this model, BHs are represented by 
collisionless sink particles initially seeded in dark matter haloes, 
which subsequently grow via gas accretion and through mergers 
with other BHs during close encounters. Every new dark matter 
halo, identified by a run-time friends-of-friends algorithm, above 
the mass threshold Mth = 10 10 Mq, is seeded with a central BH 
of initial mass 10 M®, provided the halo does not contain any 
BH yet. Each BH can then grow by local gas accretion, with a rate 
given by 



Mbh = min ( Mb , M 



(M B ,M Edd ) 



(5) 



where Mb is the accretion rate estimated with the Bondi-Hoyle- 
Littleton formula (e.g. Bondi 1952) and A/Edd is the Eddington 
accretion rate. We refer to Springel et al. (2005) for the de t ails of 
the AGN feedback implementation and to lBooth & Schavel J2009h 
for a review of the black holes accretion feedback in the context of 
cosmological simulations. An important parameter of the model is 
the radiative efficiency, e r , which gives the radiated energy L r in 
units of the energy associated to the accreted mass: 

£r = ■ (6) 
M BH C2 

Following ISpringel et ail J2005h and JFabian et al.1 fcOld) . we fix 
e r = 0.1 as a reference value, which is typical for a radiatively 



• kr37edw500: Energy-driven winds (EDW) of velocity v w — 
500 km s _1 and Kroupa IMF (the reference IMF for this work). 

• kr37mdw: Momentum-driven winds (MDW) and Kroupa 
IMF. 

• kr37agn: In this simulation the AGN feedback mechanism, 
associated to the energy released by gas accretion onto super- 
massive black holes and described in Section [3~3| is active. 

• kr37agn+edw300: In this simulation two different feedback 
models are combined: the energy-driven galactic winds of velocity 
w w = 300 km s _1 , and the AGN feedback. We will see that this 
run is completely different from the kr37agn simulation mentioned 
above. In fact, the galactic winds start to be effective at higher red- 
shift than the AGN feedback, making the haloes devoid of gas and 
consequently reducing the efficiency of the black holes accretion 
and the power of the AGN feedback. Therefore, the effect of the 
winds is more prominent than that of the AGN, with the result 
that the kr37agn+edw300 run behaves very similarly to a regular 
energy-driven winds simulation. 

• kr37nf: This simulation was run without any winds or AGN 
feedback in order to test how large the effects of the different feed- 
back prescriptions on the CIV statistics are. 

• kr37co-edw500: as we stated in Section [37X1 this simulation 
has energy-driven winds of velocity w w = 500 km s _1 , never de- 
coupled from the hydrodynamics. We will show that the coupled 
winds result in a less efficient feedback with respect to the original 
implementation. 

• ay37edw500: Energy-driven winds of velocity v v — 500 km 
s and Arimoto-Yoshii IMF. 

• sa37edw500: Energy-driven winds of velocity u w = 500 
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Figure 1. Left panel: CIV ionization fraction, as a function of the gas overdensity, at z = 3 for a gas temperature of 10 4,6 K (solid black line) and 10 50 K 
(dot-dashed green line). Right panel: CIV ionization fraction, as a function of the gas temperature, at z = 3 for an overdensity of <5 gas ~ 50 (solid black line) 
and 5 gas ~ 10 (dot-dashed green line). 



km s and Salpeter IMF (runs kr37edw500, ay37edw500 and 
sa37edw500 are identical except for the IMF). 

• kr37p400edw500: Energy-driven winds of velocity v„ — 
500 km s _1 and Kroupa IMF, with a number of particles (gas + 
dark matter) equal to 2 x 400 . This simulation is identical to the 
kr37edw500 except for the larger number of particles and was run 
in order to check for resolution effects. 



The name of each simulation reflects the choice of the IMF, the 
box size (even if the box size is equal to 37.5 h" 1 Mpc for all 
the simulati ons, we report it in order to be co nsistent with other 
works, e.g. Opp enheimer & Dav e 2006, 2008) and the feedback 
prescription adopted (only in the case of the resolution test run, 
kr37p400edw500, the number of particles is specified in the name). 
In this paper, we label as "reference runs" the first three of the list 
above: kr37edw500, kr37mdw and kr37agn. 

We use the CLOUDY code JFerlandet aJll998h to compute a- 
posteriori the CIV ioni zation fractions for e ach gas particle. As for 
the work presented in iTescari et alj ((2009), we choose the HM05 
option in CLOUDY, which consists of a UVB made by QSOs and 
galaxies with a 10% photon escape fraction and which is in agree- 
ment with other observational constraints (Bolton et al. 2005). In 
the left panel of Figure Q] we plot the CIV ionization fraction, as a 
function of the gas overdensity <5 gas + 1 = p sas /p, at z = 3 for 
a gas temperature of 10 46 K (solid black line) and 10 5 ° K (dot- 
dashed green line). In the right panel we plot the same ionization 
fraction, as a function of the gas temperature T gas , at z = 3 for an 



overdensity of S s 



50 (solid black line), this value being rep- 



resentative of the outskirts of dark matter haloes, and S s 
(dot-dashed green line). 



10 



4 STAR FORMATION RATES AND EVOLUTION OF 
THE TOTAL CIV CONTENT 

In this Section we analyse the cosmic star formation rates (SFRs) 
and the evolution of the total Carbon content in the different simu- 
lations. 

In the left panel of Figure [2] we report the total SFR of the 
simulated volume as a function of redshift. All the simulations, ex- 
cept the kr37p400edw500 run, have the same SFR down to redshift 
z ~ 9.5, when they start to differ due to their distinct feedback 
recipes. Run kr37nf (no-feedback, green dots) has the highest star 
formation rate compared to the other runs because there is no ef- 
fective mechanism able to quench the star formation. Run kr37agn 
(AGN-feedback, black dots) follows the no-feedback run kr37nf 
down to redshift 2 = 3, but then the AGN feedback becomes effec- 
tive and the SFR suddenly decreases. Runs kr37co-edw500 (cou- 
pled energy-driven winds, red dots), kr37edw500 (energy-driven 
winds, cyan dots) and kr37mdw (momentum-driven winds, blue 
dots) show more or less the same trend but for the kr37co-edw500 
the SFR is slightly higher due to the fact that the coupled winds are 
less efficient than the normal winds in making the central region of 
the haloes devoid of gas and therefore in suppressing the star for- 
mation. The kr37agn+edw300 simulation (AGN + energy-driven 
winds, yellow dots) follows the other "wind" runs (kr37co-edw500, 
kr37edw500 and kr37mdw) until redshift z = 3 but then the AGN 
feedback starts to become active and suppresses the star formation 
very efficiently. Finally, at z > 8, the kr37p400edw500 run (reso- 
lution test, grey dots), shows a higher star formation with respect 
to the other runs due to the improved resolution of this simulation 
that can resolve higher densities at earlier times. Going to lower 
redshift the SFR of the kr37p400edw500 agrees very well with that 
of the kr37edw500 run, which has the same boxsize and parameters 
but poorer resolution. This confirms the numerical convergence of 
our simulations in the regime we are interested in. We stress that 
the aim of this plot is to underline the effect of the different feed- 
back prescriptions on the SFR. For this reason observational data 
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Figure 2. Left Panel: cosmic star formation rate (SFR) for some of the hydrodynamic simulations of Table [2] Right Panel: evolution of the total Qq as a 
function of redshift for all the hydrodynamic simulations of Table[2] 



are not plotted, because our simulations can be tuned in order to fit 
the observed cosmic star formation rate at low redshift and, in this 
sense, kr37edw500 is the reference simulation: all the other runs 
in this plot differ only for the feedback prescription (all of them 
have Kroupa IMF) and were not re-tuned. We also did not plot the 
two runs with different IMFs ay37edw500 (Arimoto-Yoshii IMF) 
and sa37edw500 (Salpeter IMF), because we do not change the 
star formation efficiency as it would be required in order to match 
the observables when the number of massive stars per unit mass of 
formed stars changes: we are more interested in the chemical and 
energetic effect of the IMF. 

In the right panel of Figure [2] we show the evolution with 
redshift of the Carbon cosmological mass density fic(z)> calcu- 
lated considering the sum of the Carbon content associated to each 
particle inside the cosmological box. Energy-driven winds runs 
kr37edw500 (cyan circles— triple-dot-dashed line), ay37edw500 
(blue squares— dashed line) and sa37edw500 (red crosses— solid 
line) are identical except for the IMF, therefore the compari- 
son of these three simulations is useful to understand how Car- 
bon production is affected by the choice of the IMF. The evolu- 
tion of fic is in fact the same, but with different normalization: 
ay37edw500 has Arimoto-Yoshii IMF and produces more Car- 
bon than kr37edw500 run (Kroupa IMF) which itself produces 
more Carbon than sa37edw500 (Salpeter IMF). With the Kroupa 
IMF, the relative fraction of intermediate-mass stars, i.e. stars with 
masses 2Mq ^ m ^ 8Mq , representing i mportant Carbon pro- 
ducers, is larger than with t he Salpeter IMF dCalura & Matteuccil 
120061 ICalura & Mencill2009l) . For this reason, fi c obtained with 
the Kroupa IMF is higher than that obtained with the Salpeter IMF. 
On the other hand, the assumpti on of a top-heavy IMF, such as the 
one of lArimoto & Yoshiil dl987l) . strongly enhances the number of 
massive stars but suppresses the fraction of intermediate mass stars, 
with the overall result of a higher Carbon density at any cosmic 
epoch. Interestingly, run kr37mdw (momentum-driven winds, blue 
triangles— dotted line) follows the sa37edw500 even if the latter 



produces much less CIV as we will show in Section!?! (Figure [121 
left panel). Run kr37co-edw500 (coupled energy-driven winds, red 
squares— dot-dashed line) has the same trend as kr37edw500 but 
with a slightly higher normalization. As we stated above, this is due 
to the fact that the coupled winds are less efficient than the normal 
winds and therefore this simulation has a higher star formation that 
results in a higher Carbon abundance. This situation is even more 
extreme for the kr37nf (no-feedback, green crosses— dashed line) 
and kr37agn (AGN-feedback, black crosses— dashed line) runs, in 
fact down to redshift z ~ 3 both of these simulations do not sup- 
press effectively the star formation and for them fic is higher. Go- 
ing to lower redshift AGN feedback becomes active in quenching 
the star formation and lowers the amount of gas used by the stars 
and converts in other species (note the excess in fic with respect 
to the kr37nf run at z < 2.7). The kr37agn+edw300 run (AGN + 
energy-driven winds, yellow squares— solid line) shows the same 
fic trend of kr37edw500 even if with higher normalization at high 
redshift due to the lower strength of the winds (300 km s _1 instead 
of 500 km s _1 ). At low redshift AGN feedback becomes active, 
but, differently from the kr37agn case: i) winds have already made 
the haloes devoid of gas therefore reducing the efficiency of the 
black holes accretion and thus the power of the AGN feedback; ii) 
winds suppress star formation at high redshift lowering the Carbon 
content, so when the AGN starts to work the net result is a fur- 
ther decrease of fic . Finally, at high redshift the kr37p400edw500 
run (resolution test, black diamonds— solid line), shows a higher 
fic with respect to the kr37edw500 due to the improved resolu- 
tion of this simulation (as explained above), while at low redshift 
the Carbon content of the kr37p400edw500 approaches that of the 
kr37edw500 run. 

Comp aring the fic evolution of ou r simulations with the one 
showed in lOppenheimer & Pavel J2006I . Figure 12), we found that 
in general the normalization of our runs is higher, although a 
very precise comparison is not possible given the different setup 
of the simulations used, especially regarding the chemical evo- 
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lution model. However, for the kr37mdw run (and also for the 
sa37edw500 run) the agreement is much better, the latter simula- 
tion producing momentum-driven galactic winds and therefore be- 
ing the yardstick for the comparison with the Oppenhei mer & Pavel 
J2006h work. 



5 STATISTICS OF HI 

Before focusing on the CIV evolution, we address some flux statis- 
tics related to the properties and the evolution of neutral hydrogen 
in the IGM. For all the neutral hydrogen statistics discussed in this 
paper, we considered HI lines (while for the CIV we considered 
from time to time lines or sytems of lines). 

For each simulation performed we have extracted sev- 
eral physical quantities interpolated along lines-of-sight (LOSs) 
through the box. Given the positions, velocities, densities and tem- 
peratures of all SPH particles at a given redshift, we compute spec- 
tra along a given LO S through the box following the procedure of 
iTheunsetall dl998l) . We divide the sight line into N = 1024 bins 
of width A in distance x along the sight line. For a bin i at position 
x(i) we compute the density and the density weighted temperature 
and velocity from: 



Px 



(j) = a 3 ^X(i)W y , 



(pT)x (i) 



(Hx 0') 



a 3 ^X(i)T(0W ij: 



(V) 



(8) 



a 3 J2 X(i){ax(i) + a[x(i) ~ x(j)]}Wi h (9) 



where a is the scale factor, X(i) is the abundance of species X 
of SPH particle i, assuming ionization equilibrium, and Wy = 
mW(gij)//iy is the normalized SPH kernel. Here W is the SPH 
kernel, m is the SPH particle mass and: 

a\x(i)-x(j)\ 



<m 



his 



[Hi) + h(j)]- 



(10) 



(ii) 



with h the physical softening scale. The optical depths of the sim- 
ulated QSO spectra are drawn in redshift space taking into account 
the effect of the IGM peculiar velocities along the line-of-sight, 

^pcc, || • 

Taking the case of HI, even if the neutral hydrogen fraction 
/hi is associated to eac h gas particle and is stored in each simula- 
tion sn apshot, we follow lNagamine et alj ( 120041) . as in lTescari et al.l 
(2009), to assign a-posteriori a new mass in neutral hydrogen to gas 
particles above the density threshold for the onset of the star for- 



mation (set to pth = 0.1 cm 
hydrogen atoms) which reads: 



in terms of the number density of 



whi = fmmn (p < Pth) 
mm = /c mu (p > pth) 



(12) 
(13) 



with /hi the neutral hydrogen fraction that depends on the UVB 
used, mn the hydrogen mass of the particle (/hi and ttih are deter- 
mined self-consistently inside the code), f c the fraction of mass in 
cold clouds and pth the star formation threshold. Here f c — p c /p, 
where p c is the density of cold clouds and p the total (hot + cold) 
gas density. Individual molecular clouds cannot be resolved at the 
resolution reachable in cosmological simulations, thus p c repre- 
sents an average value com puted over small regions of the ISM 
dSpringel & Hemquistll2003l) . 



Our code follows self-consistently the evolution of Carbon 
(and of other elements, see Section |3j, then, using the CLOUDY 
code (Ferland et al. 1998), we determine the ionization fractions 
for the CIV (AA 1548.204, 1550.781 A). The simulated flux of a 
given ion transition at the redshift-space coordinate u (in km s" 1 ) 
is F(u) — exp[— t(u)] with: 



t(u) = 



i7o,i e 



iii(x) V it — x ■ 



■v^idx), b{x)\ dx, 



(14) 

where i%i is the ion number density, <to,i is the cross-section of the 
particular ion transition, H(z) is the Hubble constant at redshift z, 
x is the real-space coordinate (in km s _1 ), b — {2ksT /mic 2 ) 1 ' 2 
is the velocity dispersion in units of c, V is the Voigt profile. These 
spectra can be converted from velocity v to observed wavelength A 
using A = A (l + z)(l + v/c). 

After having extracted the spectra along random LOS through 
the cosmological box for a given simulation at a given redshift z, 
we rescaled all the HI optical depths by a constant factor Am (rang- 
ing from ~ 0.43 to ~ 0.73, depending on redshift), in such a way 
that their mean value was equal to the HI effecti ve optical depth, 
r H f , at that redshift, given by the lKim et al] J2007I) fitting formula: 



Tm = (0.0023 ± 0.0007) (1 + z) 



3.65±0.21 



(15) 



This rescaling ensures that our spectra match the observed mean 
normalized flux of the Lyman-a forest at the appropriate redshift: 
(F)m,ohs — exp(— tjh). We also rescaled all the CIV opti- 
cal depths for the same constant factor Am. Then we convolved 
the spectra with a Gaussian of 6.6 km s _1 Full Width at Half- 
Maximum (FWHM) and we add noise in order to have realistic 
spectra with signal-to-noise ratio S/N = 50, to compare with ob- 
servations. Finally, we run the VPFIT code to fit each spectra and 
recover HI and CIV lines' redshift, column density and Doppler pa- 
rameter. Among the lines fitted by VPFIT we have considered only 
those with relative errors in the Doppler parameter lower than 50%. 

5.1 The HI column density distribution function 

In this Section we compute the Lyman-a forest column density dis- 
tribution function (HI CDDF) d 2 n/dNdz, namely the number of 
absorbers with HI column density in the range [N, N + dN] and 
redshift path in the interval [z,z + dz] . The redshift path at a given 
redshift z in simulated spectra is given by: 



dz 



n sp oc(l + z) , 

c 



(16) 



where c is the speed of light (in km s ) and n spcc is the number 
of spectra taken at z for a simulation with box size in km s _1 equal 
to Av. 

In Figure [3] we show the HI CDDF at redshift z — 3 
for all the simulations of Table [2] splitted in two groups for the 
sake of clarit y. All the sim ulations fit very well the observa- 
tional data by iD'Odorico et alj (2008) down to column density 
log A^hi (cm" 2 ) = 12.2 and are all in agreement confirming the 
fact that the HI column density distribution function is quite in- 
sensitive to the different feedback prescriptions. Even though the 
galaxies drive strong winds or AGN feedbacks, there is no dis- 
cernible effect on the Lyman-a forest. This is due to the fact that 
winds or black hole ejecta expand preferentially into the lower den- 
sity regions and so keep the fil aments that produce the hydrogen 
lines intact dTheuns et al.ll2002l) . The volume filling factor of the 
winds is thus quite small and does not impact strongly on Lyman-a 
lines. In Figure|4]we show the HI CDDF at redshift z = 2.25: all 
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Figure 3. Neutral Hydrogen (HI) column density distribution function at z = 3. Left panel: part I (first set of simulations). Right panel: part II (second set of 
simulations). Data from D'Odorico et al. (2008). 
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Figure 4. As in Figure[3] but at redshift z = 2.25. 
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Figure 5. Line-widths bm probability distribution function at z = 3. Left panel: part I. Right panel: part II. In both panels, data ft om lD' Odorico et alj (2008) 
are showed by the purple solid line along with the associated Poissonian error (shaded region). 
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Figure 6. As in Figure \S\ but at redshift z = 2.25. 



the previous trends are confirmed and the agreement to the data is 
very good also for log Nm (cm -2 ) < 12. 

5.2 Probability distribution function of the HI Doppler 
parameter 

In this Section we analyse the Lyman-a Doppler parameter, &hi, 
probability distribution function. Throughout the rest of the paper, 
both the HI and the CIV Doppler parameters are assumed to be due 
only to the thermal width of the line and we neglected possible tur- 
bulent contribution. Therefore, omitting the turbulent contribution, 
the Doppler parameter basically measures the temperature of the 
gas and is defined as: 



6hi = 



2fc B r 
mm 



(17) 



where mm is the neutral hydrogen particle mass and fee is the 
Boltzmann's constant. In Figures \5\ and [6] we show the Lyman-a 
Doppler parameter distribution function at redshifts z = 3 and 
Z = 2.25, respectively, for all the simulations of Table [2] split- 
ted in two groups and compared with the data by ID' Odorico et al.l 
(2008) (purple solid line, while the shaded region represents the 
1-er Poissonian error). All runs are in good agreement with each 
other, while they are shifted towards higher 6m values with respect 
to the observational data. At redshift z — 3, the median value of 
the observational Doppler parameters distribution is 25.9 km s _1 , 
while all the simulatio ns have median around ~ 34.5 km s _1 . At 
Z — 2.25, the data bv lD'Odorico et al.l 1 12008!) have median value 
equal to 24.8 km s _1 , while the simulations have median around 
~ 34.9 km s -1 . 

Our simulations produce HI gas that is too hot compared to 
the observations. There can be many reasons for this, both physical 
and numerical. Regarding the former, it is very important to stress 
the effect of the UV background on the gas: if the UV background 
is too strong the g as will be heated too muc h. As stated in Section[3] 
we multiplied the Haardt & Madau ( 1996) heating rates by a factor 
of 3.3, introduced to fit the Schaye et al. (2000b) constraints on the 
temperature evolution of the intergalactic medium. It seems that 
this correction does not work properly with new and higher quality 
data sets. 

However, minor part of the discrepancy is due to the lines fit- 
ting procedure. When the absorption lines are weak and not well 



defined, VPFIT can add relatively broad components in order to 
minimize the \ 2 statistics. This spurious components could have 
large equivalent width and low optical depth, therefore they have 
negligible impact on the column density distribution function be- 
cause of their relative small column density, but they can bias the 
Doppler parameter probability distribution function. We will show 
later in Section[9]fhat this numerical effect is more evident for the 
CIV statistics and how it is important in explaining some features 
of the 6civ — iVciv relation. 



5.3 Column density-Doppler parameter relation 

In Figures [7] and [8] we plot the &hi — Nm relation, at redshift 
z = 3 and z = 2.25 respectively, for the reference simulations 
kr37edw500, kr37mdw and kr37agn. The contour plots are colour 
coded according to the fraction of points that fall in each bin 
with coordinate (bHi.-Vm) . In th e uppe r left panels the obser- 
vational data of iD'Odorico et alj J20081) are reported. Both at 
redshifts z — 3 and z = 2.25 all runs are in good agreement 
with each other and reproduce well the observations, even if the 
simulations present Doppler parameters slightly higher than the 
observational data (as we mentioned in the previous Section). 
The mean slope of the data envelope in the column density range 
11.5 < log Nm (c m~ 2 ) < 15 m easures the slope of the T — p 
relation for the gas dSchavduOOll) . and all the runs show roughly 
the same right slope of the data at the two redshifts. 

To sum up, the fact that for our simulations all the different 
HI statistics are in general agreement with the observational data, 
confirms that we are catching the physics of the gas traced by the 
neutral hydrogen. The discrepancies are due, mainly, to the choice 
of the UV backgroun d (assumed to be pro duced by quasars and 
galaxies as given by Haardt & Ma dau 1996, with the heating rates 
multiplied by a factor of 3.3, see Sections l3l and l5?2l and, in mi- 
nor part, to a numerical effect introduced by VPFIT that affects our 
analysis of the absorption lines in quasar spectra. We will discuss in 
more detail the impact of this effect in the next sessions regarding 
the CIV absorption lines statistics. 

The resolution-test simulation kr37p400edw500 follows the 
other runs in all the previous statistics, again confirming the nu- 
merical convergence of our simulations for the HI evolution. 
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6 THE CIV COLUMN DENSITY DISTRIBUTION 
FUNCTION 



In this Section we investigate the CIV column density distribu- 
tion function (CIV CDDF, already defined in Section [5~lt . plotted 
at redshift z = 3, 2.2 5 and 1.5, respectively and along with the 
ID'Odorico et all d2010l) data, in Figures [9] \W\ and Qj] We consid- 
ered for this analysis systems of lines as defined in Section[2] 

At z - 3 (Figure[9}, runs (left panel) kr37agn+edw300 (AGN 
+ energy-driven winds), kr37co-edw500 (coupled energy-driven 
winds), kr37mdw (momentum-driven winds) and (right panel) 
kr37edw500 (energy-driven winds, Kroupa IMF), ay37edw500 
(energy-driven winds, Arimoto-Yoshii IMF), sa37edw500 (energy- 
driven winds, Salpeter IMF), are in agreement with each 
other and reproduce well the observed distribution down to 
log iVcrv (cm -2 ) ~ 12.7. In contrast, for the simulations kr37nf 
and kr37agn the agreement is not good. Simulation kr37nf (Fig- 
ure [9] right panel, green dashed line) was run, as a test, without 
any winds or AGN feedback: so there is no effective mechanism 
to spread the enriched gas outside the haloes. The result is that 
all the gas is trapped inside high density haloes and most of the 
Carbon is locked back into stars. In the run kr37agn (Figure [9] left 
panel, black dashed line) the feedback mechanism is the AGN feed- 
back, associated to the energy released by gas accretion onto super- 
massive black holes. Comparing the resulting distribution with that 
of kr37nf, we note that they are very similar. The reason is that at 
redshift z — 3 the AGN feedback is not active yet, but it is just start- 
ing to work, so at high redshift this simulation virtually behaves like 
the no-feedback one. Interestingly, run kr37agn+edw300 (AGN + 
energy-driven winds) follows all the other "wind" runs even if in 
this case we have the combined effect of winds and AGN feed- 
backs: this is due to the fact that winds start to be effective at higher 
redshift than the AGN feedback, making the haloes devoid of gas 
and in this way reducing the efficiency of the black holes accretion 
and the power of the AGN feedback. 

At redshift z — 2.25 (Figure I10t . the agreement with ob- 
servational data is very good down to log /Vciv (cm -2 ) = 12. 
The kr37nf run (no-feedback, right panel dashed green line) suf- 
fers of the same problems described before, while for the kr37agn 
run (AGN feedback, left panel dashed black line) the discrepan- 
cies with the observational data decreased: a first hint that going 
to lower redshift AGN feedback becomes more and more efficient 
and the behaviour of this simulation starts to resemble that of the 
"wind" runs. 

At redshift z = 1.5 (Figure [TT) all the simulations produce 
too many systems with low column density and systematically fall 
slightly short of reproducing the distribution at log iVcrv (cm -2 ) > 
14.5. Now the kr37agn run (AGN feedback, left panel dashed black 
line) closely follows all the other runs, confirming that at this low 
redshift AGN feedback is active and powerful, while kr37nf run (no 
feedback, right panel dashed green line) still fails in reproducing 
the correct trend. For this statistics, at the three redshift considered, 
the resolution-test simulation kr37p400edw500 follows the refer- 
ence runs. 

Note that at redshift z — 2.25 we fit very well the obser- 
vational data in the low column density tail of the distribution, 
log iVcrv (cm -2 ) < 12.2, while at z = 3 and 2 = 1.5 the agree- 
ment in this column density range gets worse. 



7 THE EVOLUTION OF THE CIV COSMOLOGICAL 
MASS DENSITY 

In the left panel of Figure [T2] we show Qciv(z), calculated 
considering the sum of the CIV content associated to each 
particle inside the cosmological box. Energy-driven winds runs 
kr37edw500 (cyan circles— triple dot-dashed line), ay37edw500 
(blue squares— dashed line) and sa37edw500 (red crosses— solid 
line) behave exactly as we described in Section [4] they differ only 
for the IMF and the ay37edw500 produces much more CIV than 
the other two, having Arimoto-Yoshii IMF. Runs kr37co-edw500 
(coupled energy-driven winds, red squares— dot-dashed line) and 
kr37mdw (momentum-driven winds, blue triangles— dotted line) 
follow almost the same trend as kr37edw500, which has the same 
IMF of the two runs mentioned above. For the kr37nf simulation 
(no-feedback, green crosses— dashed line), ficiv starts to decrease 
at redshift z ~ 2.7. It is interesting to note that for the simulation 
kr37agn (AGN-feedback, black crosses— dashed line) ficiv sud- 
denly increases after a small decrement at redshift z ~ 2.8 (down 
to this redshift the kr37agn run follows the no-feedback run). This 
is due to the suppression of star formation produced by the power- 
ful AGN feedback. Then, moving to lower redshift the CIV cosmo- 
logical density starts to approach all the other runs. The same effect 
is visible in the kr37agn+edw300 run (AGN + energy-driven winds, 
yellow squares— solid line). In this case the effect is weaker be- 
cause, as we explained in Section[6] winds make the haloes devoid 
of gas therefore reducing the efficiency of the black holes accretion 
and thus the power of the AGN feedback. As for the total S7c plot 
in the right panel of Figure [2] at high redshift, the resolution-test 
simulation kr37p400edw500, shows an higher £7civ with respect 
to the kr37edw500, while the two runs are in agreement at low red- 
shift. 

In the right panel of Figure Q/2] we show the evolution with 
redshift of the CIV cosmological mass density as a fraction of the 
critical density today, calculated considering the column density 
distribution function of the CIV absorption lines: 



^civ(z) = 



H mciv 

C Pcrit 



Nf(N) dN, 



(18) 



where Ho = 100 h km s 1 Mpc : is the Hubble constant 
(h — 0.73), mciv is the mass of a CIV ion, c is the speed of light, 
p cxix = 1.88xl0 -29 h 2 gem -3 and f(N) is the CIV CDDF. Since 
f(N) cannot be recovered correctly for all the column densities due 
to incompleteness and poor statistics, the integral in the previous 
equation can be approximated by a sum dStorrie-Lombardi et al.1 
1 1996b : 



iiciv(z) — 



HomcwEiNiiCIV) 



C Pcrit 



AX 



-1/2 



(19) 



dz, the red- 



Wifh AX = J(l + zf [fi rn(l + zf + Q a] 

shift absorption path. In the right panel of Figure [T2] the overplot- 
ted black diamonds show the observational data of D'Odorico et al. 
(2010). The value of Slciv significantly depends on the column 
density range over which the sum or the integration is carried 
out and, as a consequence, on the resolution and signal-to-noise 
ratio of the ava ilable spectra. In order to address this aspect, 
ID'Odorico et alj (2010) have computed three sets of values to be 
compared consistently with different data in the literature. Here 
we compared Ociv obtained from the simulated CIV systems (as 
defined in Section |2) lying in the column density range 12 ^ 
log iVciv (cm -2 ) ^ 15, with the ^civ co mputed from the CIV 
systems of the total ID'Odorico et all (2010) sample, in the same 
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Figure 9. CIV column density distribution function at z = 3. Left panel: part I. Right panel: part II. Data from lD'Odorico et alji201dh . 
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Figure 10. As in Figure|5] but at redshift z = 2.25. 
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Figure 11. As in Figure [9] but at redshift z = 1.5 
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range of column densities. From z — 3.5 to z = 2.5 all the simula- 
tions, except for the kr37nf (no-feedback) and the kr37agn (AGN- 
feedback), roughly reproduce the observational data and show an 
increasing or constant trend for the ficiv • Down to redshift z — 3 
the AGN feedback of run kr37agn is not active and the CIV content 
is nearly equal to the one of the no-feedback simulation kr37nf. 
At z = 3 the AGN feedback starts to work and the kr37agn run 
suddenly reaches all the other "wind" simulations. 

Apart from the no-feedback simulation kr37nf (green crosses 
and dashed line), all the other runs still reproduce the observational 
data around redshift z = 2.25 ± 0.25. If we consider lower red- 
shift, all the simulations show a decreasing trend at variance with 
the data. Different is the case of the kr37agn run: from redshift 
z = 3.0 to 1.8, the AGN feedback suppresses efficiently the star 
formation, therefore the gas is no longer reprocessed from the stars 
and the CIV is not converted in other ions. At the same time, the 
AGN feedback has not yet made the haloes devoid of gas, with the 
result that firjiv continues to increase. However, since the black 
holes accretion feedback is extremely powerful, moving to lower 
redshift a considerable amount of gas is expelled from the haloes 
and heated. As a consequence Q.cw starts to decrease, even if at a 
later time and to a smaller degree than the "wind" runs. 



8 PROBABILITY DISTRIBUTION FUNCTION OF THE 
CIV DOPPLER PARAMETER 

In this Section we focus on the CIV Doppler parameter probability 
distribution function, plotted at redshift z = 3, 2.25 and 1.5 in 
Figure s l 1311 14l and| 151 respectively. The data from lD'Odorico et al.1 
(2010) are also overplotted in these figures with a purple solid line 
along with the associated Poissonian error (shaded region). For this 
analysis we considered CIV lines. 

At redshift z = 3 (Figure [T3l all the runs are quite in agree- 
ment with the data even if there are some small discrepancies. All 



the runs tend to underestimate the observed distribution in the first 
bin, bciv < 4 km s _1 , and slightly overestimate it for &civ > 28 
km s _1 . In the intermediate range 4 < &civ < 28 km s _1 , the 
simulations behave differently. The kr37mdw, kr37agn+edw300, 
kr37co-edw500 (left panel) and sa37edw500 (right panel) runs 
fit very well the observational data, while runs kr37edw500 and 
ay37edw500 (right panel) overproduce the observed distribution in 
the bin centered at &civ = 10 km s _1 . Again the runs kr37nf and 
kr37agn show a different trend with respect to the others, for ex- 
actly the same reasons mentioned in the previous sections. Sim- 
ulations kr37nf (no-feedback, right panel green dashed line) and 
kr37agn (AGN-feedback, left panel black dashed line) distributions 
are shifted towards higher bciv than the observed one. At this red- 
shift AGN feedback is not active so, as the no-feedback case of 
kr37nf, the gas remains trapped inside haloes at high density and 
temperature. Run kr37p400edw500 show a marked excess of low 
CIV Doppler parameters due to its improved resolution that can re- 
solve higher densities at earlier times and produces a lot of small 
metal enriched clumps in the IGM. At z — 3, the median value 
of the observational Doppler parameters distribution is 9.88 km 
s _1 , while the simulations have median ranging from 8.49 km s _1 
(kr37p400edw500) to 13.03 km s~ x (kr37nf). 

The situation is different at redshift z = 2.25 (Figure [l4"t 
and z — 1.5 (Figure [T5l>. Moving to lower redshift, all the sim- 
ulations show an excess of low Doppler parameters (6civ < 8 
km s _1 ), and also they underproduce the observed distribution in 
the range 8 < &civ < 20 km s~\ Moreover, while at low red- 
shift the kr37nf simulation (no-feedback, right panels green dashed 
lines) is still shifted towards higher &civ with respect to the ob- 
served distribution, the kr37agn run (AGN-feedback, left panels 
black dashed lines) at redshift z — 2.25 fits better the data and 
it finally approaches all the other runs at z — 1.5, further confirm- 
ing the increasing efficiency of the A GN feedback at low re dshift. 
At z = 2.25 and z = 1.5, the data bv lD'Odorico et alj j2oTob have 
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Figure 14. As in Figure [Til but at redshift z = 2.25. 
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Figure 15. As in Figure [T3l but at redshift 2 = 1.5. 
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median value respectively equal to 9.30 and 9.25 km s~ , while the 
simulations have median around, respectively, 8.55 and 7.40 km 
s _1 (12.47 km s _1 in the case of the no-feedback run kr37nf). 

We stress that the small discrepancy between observed and 
simulated distributions for low Doppler parameters (&civ < 8 km 
s _1 ) appears in a range that is influenced by the adopted resolution 
(remember that we have convolved the simulated spectra with a 
Gaussian of 6.6 km s _1 FWHM) and a sightly different value than 
the one used could alleviate this discrepancy 

Compared to that of HI, the CIV Doppler parameter distribu- 
tion function is less affected by the choice of the UV background 
(and in particular of the factor of 3.3 discussed in Sections [3] and 
15.2b . This is due to the fact that CIV absorption is related to higher 
density regions (mainly the outskirts of galactic haloes) than the HI 
absorption. In these regions, the impact of the IGM local temper- 
ature and density on the ionization state of Carbon is considerably 
more important than the UVB contribution. 



9 THE CIV COLUMN DENSITY-DOPPLER 
PARAMETER RELATION 

In Figures[T6][T7]and[T8]we plot the feciv — iVcrv relation, at red- 
shift z = 3, z = 2.25 and z — 1.5 respectively, for the reference 
simulations kr37edw500, kr37mdw and kr37agn. The contour plots 
are colour coded according to the fraction of points that fall in each 
bin with coordinate (bciv^-Vc iv). In the upper left pan el of the fig- 
ures the observational data of iD'Odorico et al] 1 120 101) are plotted. 
As in the previous Section, we considered here CIV lines. 

At redshift z — 3 (Figure [Tol l, the data distributions of runs 



kr37edw500 (energy-driven winds) and kr37mdw (momentum- 
driven winds) are in agreement with each other but are quite differ- 
ent from the observed one. In particular the simulated distributions 
are shifted towards higher CIV column densities and are concen- 
trated around a narrower region than the observational data. The 
low column density tail of the distributions is completely missing, 
while the high column density tail shows a sort of correlation be- 
tween the CIV Doppler parameters and Nctv. We interpret this 
latter as a T — p relation: lines with higher column densities (i.e. 
associated to denser regions) have higher Doppler parameters (i.e. 
temperature). The fact that this relation appears only at high column 
densities is not a surprise: for high values of iVciv. absorption lines 
are very strong and well resolved, so they are good tracers of the 
physical state of the intergalactic gas and also they are better fit by 
VPFIT. The case of run kr37agn is different: in fact the T — p re- 
lation for large Nciv values is missing, while a clump of systems 
with a large spread in feciv values at low column density is clearly 
visible. This distribution results because, at redshift z = 3, in the 
kr37agn run the gas resides inside the high density high tempera- 
ture cores of the haloes where it is reprocessed by the stars, and the 
amount of diffused CIV is rather low (see the right panel of Fig- 
ure!12t. Therefore, since the CIV absorption lines are weak and not 
well defined, the quality of the fit made by VPFIT is rather poor. As 
we reported in Section [5T2l in such a case VPFIT tends to add broad 
components in order to minimize the \ 2 statistics and this results in 
the feature of spurious systems at low column density shown in the 
lower right panel of Figure[l6] While the effect is less evident (but 
present) for the other runs or in the case of the bui — Nm relation, it 
is now much more prominent. Of course, observational data do not 



© 0000 RAS, MNRAS 000, 000-000 



Cosmic evolution of the CIV 17 



0.00 0.01 0.02 0.03 0.04 0.05 



0.01 0.02 0.03 0.04 0.05 



22.5 



& 15 



7.5 



D'Odorico el at. (2010) 





0.00 0.01 0.02 0.03 0.04 0.05 



0.01 0.0. 



0.04 0.05 



22.5 



rfi 15 



7.5 



11 




12 13 

log 10 N CIV 



14 11 




12 13 

log io Nav 



Figure 17. As in Figure [T6l but at redshift z = 2.25. 
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Figure 18. As in Figure [T6l but at redshift 2 = 1.5. 
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present such a feature, because in the "by eye" procedure of fitting 
these spurious lines are removed. 

At redshift z — 2.25 (Figure I17t and z = 1.5 (Figure 
|18| ) the kr37agn run approaches the other two simulations. The 
&civ — -/Vciv correlation at high column densities is now present 
for all the runs but the simulated distributions are still different 
from the observed one. The clump of low density systems disap- 
pears almost completely at z = 2.25, but is visible also for the runs 
kr37edw500 and kr37mdw at redshift z = 1.5. There is a simple 
reason for this: at low redshift the amount of CIV, traced by ran- 
dom line-of-sight along the box, decreases (right panel of Figure 
I12t , while the number of low column density systems significantly 
increases (Figure fTTTl and typically these lines have a large spread 
in their Doppler parameter value. 

The difference between simulated and observed distributions 
could also be due to the fact that simulated haloes, and in particular 
their outskirts, can have different properties from that of the ob- 
served ones: for example this is also seen from the Sill statistics in 
DLA systems that appear to be not in agreement with observations. 
Although in the present paper we explore the effect of feedback, 
there might be other physical effects that can contribute as well in 
shaping the properties of the &civ — A^civ distribution. 



10 Afciv - iVm CORRELATED ABSORPTION 

The last part of this work is dedicated to the HI-CIV correlated 
absorption. We define as correlated absorption HI-CIV systems in 
which CIV and HI are physically dependent or, in other words, 
related to the same absorptive structure. In real spectra, however, 
two absorbers each containing CIV could be very close in redshift 
space but far in real space due to the bulk motions and the peculiar 
velocity. For this reason there are often serious problem on how 
to associate each CIV components to physically corresponding HI 



2010), w e comp are with the fitting formulae obtained from the 



compo nents. In this Section, besides the data of ID'Odorico et al.l 

ive comp a 
iKim et al.l {20 10) data. In oder to do that, only for this analysis we 
changed the velocity range, du m in, used so far to define systems of 
lines. Following Kim et alj (2010) and differently from what was 
done in Sections [6] and [7] where we set df m i n = 50 km s -1 , we 
use here the velocity interval [-250 km s -1 , 250 km s -1 ] both for 
the simulated and the observed CIV and HI systems. 

In Figures [19] and [20] we s how our findings for s ome of the 
simulations in Tableland for the lD'Odorico et al. ( 2010) da ta. The 
red-dashed and yellow-solid curves represent the lKim et al.l q2010h 
fit, obtained by using the function: 



log New = 



d 



log N m + C 2 



+ C*3. 



(20) 



The parameters for the red-dashed curves are: G\ — [—4.43 ± 
1.30,-7.57 ± 1.80], C 2 = [-13.04 ± 0.44,-12.45 ± 0.45], 
C-i = [15.22 ± 0.17, 15.65 ± 0.19] at z = [2.25 3.0] These red- 
dashed curves represent a fit to all the lKim et alj (2010) data avail- 
able at the two different redshifts. The parameters for the yellow- 
solid curves are: Ci = [-2.01 ± 0.35, -2.51 ± 0.47], C 2 = 
[-14.22±0.14, -14.08±0.20], C 3 = [14.96±0.05, 15.09±0.07] 
at z — [2.25, 3.0]. These curves represent an improved fit made to 
reproduce better the clump of systems at lo w Nm and low Ngiv- 
At redshift z = 3 (Figure QUl for the lD'Odorico et all j201Ch 
observational data and for the kr37agn+edw300, kr37edw500 and 
kr37mdw runs the bulk of systems is located at log AThi (cm -2 ) < 
16 and it follows fhe iKim et al. (20 10) fitting functions, even if with 



a slightly steeper trend, while at log AThi (cm - ) > 16 the spread 
in the simulated data increases (also the data of IKim et alj|2010l 
present a sort of bimodality with high HI column density distribu- 
tion more spread than the tight low HI column density one). Runs 
kr37nf (no-feedback) and kr37agn (AGN-feedback) at this redshift 
behave similarly: they do not show the log Am (cm -2 ) < 16 
CIV-HI correlation above, and the logAm (cm -2 ) > 16 sys- 
tems are rare. All the simulations show a strip of values with 
log A/hi (cm -2 ) < 15 and log Nciv (cm -2 ) ~ 12. This is due to 
the VPFIT spurious line fitting effect described in the previous Sec- 
tion. In fact, it is particular evident for the runs kr37nf and kr37agn, 
because in these two cases the amount of the diffuse CIV is lower 
and the CIV absorption lines are weak and not well defined. 

At redshift z = 2.25 (Figure I20t the above trends are in 
place even if all the simulated distributions are now shifted towards 
slightly lower HI column densities with respect to the fitting func- 
1 2010). Otherwise, the number of components 



tions of 



Kim ct al. 



of the ID'Odorico et alj ( 1201 0) data increases at this redshift and 
the bulk of systems at log N m (cm -2 ) < 16 is well fitted by the 
functions o flKim et alj 1201U). Again, moving to lower redshift, the 
simulation with AGN feedback, kr37agn, starts to follow the other 
wind runs. 

To summarize, two effects are visible from the plots: (a) for 
a given AThi, the scatter in the corresponding New is broader, 
and the Aciv values typically larger, than what is inferred from 
observations; (6) the s imulated d i stribu tions of Nm are slightly 
lower compared to the IKim et alj ( 120 101) fitting functions and the 
ID'Odorico et"aU d2010l) data at redshift z = 2.25. 



11 CONCLUSIONS AND DISCUSSION 

In this paper we investigated the global properties of triply ionized 
Carbon (CIV) as an IGM tracer at z ^ 1.5. From the numerical 
point of view we presented the results of a new set of hydrodynamic 
simulations that incorporate feedback either in the form of galactic 
winds or in the form of energy from black holes accretion (AGN 
feedback). From the observational point of view we relied on recent 
high-resolution data sets obtained with the UVES spectrograph at 
VLT and the HIRES spectrograph at Keck. The main results can be 
summarized as follows: 

• The statistics of HI are weakly affected by the feedback pre- 
scription implemented: both the column density distribution func- 
tion and the Doppler width distributions do not change significantly 
when considering the different simulations. We regard this as a ge- 
ometrical effect: winds and AGN feedback are stronger at the inter- 
section of IGM filaments and these sites have usually small filling 
factor. Only the high column density tail shows small differences 
between the runs considered (Section[5j- The fact that overall our 
simulations reproduce well all the HI statistics confirms that we are 
catching the physics of the gas traced by the neutral hydrogen. The 
remaining small discrepancies with th e observational data co uld be 
due to the UV background we used JHaardt & Madavj|l99q . pro- 
duced by quasars and galaxies, with the heating rates multiplied 
by a factor of 3.3), that heats too much the gas. In minor part, the 
discrepancies are also due to a numerical effect introduced by the 
line-fitting software VPFIT. 

• Compared to that of HI, the statistics of CIV is strongly af- 
fected by feedback. If we consider systems of lines as defined in 
Section[2] at redshift z = 3, 2.25 and 1.5, all the simulations except 
for the no-feedback run kr37nf and the AGN-feedback run kr37agn 
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Figure 19. Column densities Nay ~~ ^HI correlated absorption for s ystems at 2 = 3. Upper left panel: observational data from D'Odorico et al. (2010). 
The red-dashed and yellow-solid curves represent the fitting functions of Kim et al. 12010). 



(this latter only at high redshift) agree with the observed CIV col- 
umn density distribution function. Galactic winds feedback starts 
to be active at high redshift, but moving to lower redshift, also the 
AGN feedback becomes effective and the kr37agn simulation fol- 
lows all the other wind runs (see Figuresl9lll0land[TTV 

• In the second part of Section [7] we presented the evolution 
with redshift of the CIV cosmological mass density, calculated con- 
sidering the column density distribution function of the CIV ab- 
sorption lines. In the redshift range 2 = 2.5 — 3.5, our simula- 
tions, except for the no-feedback run kr37nf and the AGN-feedback 



run kr37agn, reproduce the observed ficrv(z) evolution, even if 
with a slightly overestimation, and they perfectly reproduce the 
observational data around redshift 2 = 2.25 ± 0.25. At lower 
redshift we found a decreasing trend at variance with the increas- 
ing trend shown by the observational data (right panel of Figure 
I12t . In the AGN case the trend is nearly constant at low redshift 
and in better agreement with data. We think that a better calibra- 
tion of the coupled AGN + galactic winds feedback model of run 
kr37agn+edw300 could help to improve this result. 

• At all the redshift considered, the CIV Doppler parameter dis- 
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Figure 20. As in Figure[l9] but at redshift z = 2.25. 



tribution is in good agreement with the observational data. The no- 
feedback and AGN-feedback simulations (the latter only at high 
redshift) result in distributions shifted towards higher bciv than 
the observed ones: this is due to the fact that these runs produce gas 
around galaxies that is too hot than what was inferred from obser- 
vations (Section[8}. 

• The overall &civ — Nciv data distribution of our reference 
simulations is different from the observed one. The high column 
density tail of all the runs considered (but for the AGN-feedback 
run this is true only at low redshift) shows a sort of correlation 
between the CIV Doppler parameters and the column density, in 



which the Doppler parameters increase as iVcrv increases: the sign 
of a well defined temperature-density relation for the high density 
gas. Especially at low redshift, when the amount of CIV traced by 
random line-of-sight along the box decreases, a clump of systems 
with a huge scatter in the &civ values at low CIV column density is 
present. This effect should be regarded as of numerical origin and it 
is introduced by VPFIT when fitting regions close to the continuum 
in the mock data set (Section [9}. A more careful comparison be- 
tween data and simulations is needed in order to address this effect 
at a more quantitative level. 

• In Section [10] we explored the correlated CIV-HI absorption, 
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considering systems of lines in which CIV and HI are physically de- 
pendent. At the two redshift considered, z = 3 an d z = 2.25, both 
our re ference runs and the observational data of D'Odoric o et al.l 
(2010) roughly follow the fit proposed bv lKimet alj J201CJ) . even 
if the simulations show a slightly steeper trend in the range 
logiV H i(cm- 2 ) < 16. 

• Even if we are not able to fit all the statistics at once if we con- 
sider the neutral hydrogen HI and the triply ionized Carbon CIV, 
momentum-driven winds simulation reproduce best all the different 
quantities explored in this paper (as found for Damped Lyman-a 
systems by Tescari et al. 2009, and for IGM metal lines by Oppen- 
heimer & Dave 2006, 2008). 

• To sum up, feedback appears to be a crucial physical ingre- 
dient in order to reproduce statistics of metal absorption lines. In 
particular, the effect of galactic winds triggered by supernovae is 
important at high redshift, while the effect of black holes accretion 
feedback starts to become prominent at z < 3. However, because 
the AGN feedback is very strong, it will be crucial to improve the 
coupled AGN + galactic winds model, where the effect of the winds 
at high redshift reduces the strength of the black holes accretion at 
low redshift. 

In our analysis, we found two different effects biasing the 
statistics related to the HI and CIV absorption lines: the choice of 
the UV background and the VPFIT lines fitting procedure. In the 
future we will try to improve these statistics by using a more physi- 
cally motivated UV background that takes into account Hell reion- 
ization at z ~ 3. Furthermore, we will try to improve our automatic 
line-fitting analysis by removing numerical artifacts introduced by 
VPFIT, since the automatic fitting procedure is not likely to be as 
accurate as a "by eye" fitting made by an observer, and introduce 
spurious lines especially in the continuum region. Our results can 
also be tested by using other statistics based on pixel-optical depths 
techniques. 

A second possible improvement of the results presented in 
this paper is the refinement of the feedback models and the in- 
clusion of new physical processes like turbulence, metal diffusion 
and radiative transfer, that can help in polluting with metals the 
low-density IGM. So far we explored feedback prescriptions in the 
form of energy-driven galactic winds (EDW), momentum-driven 
galactic winds (MDW) and black holes accretion (AGN feedback). 
Exploring the combined effect of these models seems to be promis- 
ing: for e xample the "energy coup led momentum driven winds" 
model bv lChoi & Nagamind J2010]) succeeds in enriching the gas 
(like the EDW model) without heating too much the IGM (like the 
MDW model). The modelization inside the hydrodynamic simula- 
tions of the radiative transfer is fundamental and could also impact 
on the metal mixing. Finally, it will be crucial to incorporate in the 
GADGET-2 code the small-scale turbulence and its impact on the 
metal diffusion at large scales. With physical mixing, fluid elements 
on a fixed (resolved) physical scale do exchange energy/entropy 
due to unresolved (turbulent) motions: diffu sion allows so me ejecta 
gas to mix while exiting the galaxies iShenet alj|2010f) . All the 
previous effects, once included in simulations, could help in im- 
proving our understanding of the chemical and physical evolution 
of the IGM and provide a more comprehensive framework of the 
high-redshift galaxy/IGM interplay. 
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